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Abstract 

The effects of local electronic interactions and finite temperatures upon the transmission across 
the CU 4 C 0 CU 4 metallic heterostructure are studied in a combined density functional and dynamical 
mean held theory. It is shown that, as the electronic correlations are taken into account via a local 
but dynamic self-energy, the total transmission at the Fermi level gets reduced (predominantly in 
the minority spin channel), whereby the spin polarization of the transmission increases. The latter 
is due to a more signihcant d-electrons contribution, as compared to the non-correlated case in 
which the transport is dominated by s and p electrons. 
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I. INTRODUCTION 


The design of multi-layered heterostructures composed of alternating magnetic and non¬ 
magnetic metals offers large flexibility in tailoring spin-sensitive electron transport properties 
of devices in which the current flow is perpendicular to the planes. In the framework of ballis¬ 
tic transport, the spin-polarized conductance and the giant magnetoresistance effect (GMR) 
depends on^h^ mismatch between the electronic bands of the concerned metals near the 

arization of current and hence the GMR, 


Fermi level [l, 2|. In order to maximize the spin po’ 
heterostructures including half-metallic materials 


^ seem to be the materials of choice. 


In practice, however, the spin-polarization is never complete due to the presence of defects. 


Q. 


and/or due to intrinsic limitations caused by spin-contamination and spin-orbit coupling 
Owing to the technological relevance, considerable progress has been achieved in the com¬ 
putational description of multilayered heterostructures. In particular, the ballistic transport 
properties have been addressed by considering the Landauer-Biittiker formalism jb-O], where 
the conductance is determined by the electron transmission probability through the device 
region, which is placed between two semi-infinite electrodes. The transmission probability 
can be then computed with different electronic structure approaches, such as the tight- 
binding 10l-ll3| or the first-principles density functional theory (DFT) ones 14j-ll6|. Various 
implementations exist, based on transfer matrix [u, 18|, layer-Korringa-Kohn-Rostoker 
(KKR) [l9, 201, or non-equilibrium Green’s function (NEGF) 211 techniques. 

In the context of first-principles calculations, it is known that, for systems with moderate 
to strong electron correlations, the electronic structure, calculated with the “conventional” 
DFT local density approximation (LDA) or its generalized gradient approximation (GGA) 
extension, is not accurate enough to account for the observed spectroscopic behavior. A more 
adequate description is provided within the dynamical mean field theory (DMFT) 221-124 1 


built on the LDA framework 


25 


27l |. Since many interesting magnetic materials fall into 


this category, the prediction of their electron transport properties, obtained by combining 


the Landauer-Biittiker formalism with DFT 


21 


28 


29|, is expected to equally suffer from 


an insufficient treatement of correlation effects, which would be captured by adding DMFT. 
However, the full incorporation of correlation effects at the DMFT level into the transport 
calculations is not straightforward, not only because of technical reasons, such as large sys¬ 
tem size and lack of corresponding algorithms, but also because of conceptual difficulties in 
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the development of many-body solvers in the out-of-equilibrium regime 3^. Attempts to 
close this gap include the use of combined techniques, in which equilibrium-DMFT calcula¬ 
tions are performed in order to obtain the Landauer conductance of atomic contacts made 
of transition metals j^. 32|. 

In this work, we investigate the linear-response transport through a prototypical Cu- 
Co-Cu heterostructure by accounting for strong electron correlation effects in the electronic 
structure of the Co monolayer. The attention to this system is drawn by a signihcant density 
of states that develops in the Co layer, in the vicinity of the Fermi level, in one spin channel 
only (the minority-spin one), whereas the Cu layers contribute with states at higher binding 
energies. We use a “two-step” approach, in which the Landauer transmission probability is 


calculated within the smeagol NEGF based electron transport code 


Hamiltonian is obtained from DFT 


Hfl 


29j whereby the 


33| . The many-body corrections to the Green’s func 


tion 


are evaluated using DMFT in an exact muffin-tin orbitals (EMTO)-based package 


34 


36|, 


which uses a screened KKR approach [37[. These corrections are then passed to smeagol 
for the calculation of the transport properties. 

The article is organized as follows. We start with a general description of the transport 
problem in the presence of electronic correlations fSec. Ill A|) . Then the computational details 
are outlined in Sec. Ill Bl and the geometry of the system considered in our simulations is 
descried in Sec. Ill Cl Finally, Sec. IIIII presents the main results, and Sec. HYl summarizes and 
concludes. The appendices deal with technical implementation and various tests, notably 
App. |B] explains porting the results from EMTO into smeagol, and discusses a model 
two-orbital system of cubic symmetry. 


II. METHODS 


A. Transport properties in the presence of electronic correlations 


The electronic transport through a device can be addressed using the Kubo approach, 
where the central quantity is the conductivity, and the electrical current is the result of 


the linear response of the system to an applied electric held 38|. Alternatively, in the 


Landauer-Biittiker formulation jb-^, the current how through a device is considered as a 
transmission process across a hnite-sized scattering region placed between two semi-inhnite 
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leads, connected at infinity to charge reservoirs. The quantity of interest is the conductance, 
which, within linear response, is given by: 

e 


g = 


.2 1 


h 


BZ 


V [ dk||T.(k||,E^), 


( 1 ) 


where —e is the electron charge, h the Planck constant, and e^/h half the quantum of 
conductance. To-(k||, Ep) is the spin-dependent transmission probability from one lead to the 
other for electrons at the Fermi energy and with the transverse wave-vector ky perpendicular 
to the current flow (here we assume that the two spin-components do not mix). The integral 
over k|| goes over the Brillouin zone (BZ) perpendicular to the transport direction, and f^Bz is 
the area of the BZ. In the case when the interaction between electrons involved in transport 
is completely ne glec ted, the transmission for a given energy, E, of the incident electrons can 


be evaluated as 


39| 


T,(k||, B) = [ri(k||, B)G'’t(k||, B)rj(k||, B)G'’(k||,B)] , 


( 2 ) 


where G'^(k||,i?) is the retarded Green’s function of the scattering region coupled to the 
leads, 

G-'(k||,B)= [€+S(k||)-H''(k||)-SI(k||,B)-S;,(k||,i5)]-‘. (3) 

All terms presented are matrices [G'’'(k||, E)]^^, labelled by the global indices fi, v which run 
through the basis functions at all atomic positions in the scattering region. S(k||) represents 
the orbital overlap matrix, and the energy shift into the complex plane, e"*" = hm 5 ^o+ {E+i6), 
has been introduced to respect causality. H°'(k||) is the Hamiltonian of the scattering 
region for spin a; the right and left self-energies S^(k||,F^) and S2(k||,ii^) describe the 
energy-, momentum- and spin-dependent hybridization of the scattering region with the 
left and right leads, respectively 29|]. Therefore, G°'(k||,i?) is formally the retarded Green 
function associated to the effective, non-hermitian Hamiltonian Hgg(k||,ii^) = H°'(k||) — 
SI(k||.i!)- S^(k||,i?), in which the self-energies act as external energy-, momentum- and 
spin-dependent potentials. In Eq. ([2]), r2(^)(k||,E) = i [S2(^)(k||, E) — S^|^^(k||, E)] is 
the so-called left (right) broadening matrix that accounts for the hybridization-induced 
broadening of the single-particle energy levels of the scattering region. Importantly, for 
non-interactin g el ectrons, it has been proved that the Landauer and the Kubo approaches 
are equivalent ^ , so that the linear-response transport properties of a system can be com¬ 
puted with either formalism. During the last few years, the Landauer approach has been 
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systematically applied in conjunction with DFT in order to perform calculations of the con¬ 


ductance of different classes of real nano-devices 


4l| . In this combination the DFT provides a 


single-particle theory in which the Kohn-Sham eigenstates are interpreted as single-particle 
excitations. Although this approach is only valid approximatively, DFT-based transport 
studies have provided insightful results concerning the role of the band-structure in the 


electron transport process through layered heterostructures 


flS, 


42 


44 |. 


With the effect of the electron-electron interaction beyond the DFT explicitly considered, 
the retarded Green’s function of Eq. ([3]) is replaced by the following one, carrying the 
subscript “MB” for “many-body”: 


Gifl(k||,B) = [£+S(k||) - H”(k||) - SJ(k||,B) - EJ(k||,B) - SX,B(k||,B)] 


-1 


( 4 ) 


Here, S^^(k||,E') is the many-body self-energy dehned through the Dyson equation 




This accounts for all electron-electron 


^Ms(k||,-E) — G°'(k||,E) GXfB(k||,E) 
interaction effects neglected in G°^(k||,E). The self-energy acts as a spin-, momentum- 
and energy-dependent potential, whose imaginary part produces a broadening of the single¬ 
particle states due to hnite electron-electron lifetime. In this work, the many-body self¬ 
energy is computed at the DMFT level, meaning that S^^(k||, E) is approximated by a k- 
independent quantity i.e., a spatially local but energy-dependent potential. Then, 

as suggested by Jacob et al. {m, 3; conductance and the transmission probability 
are obtained within the Landauer approach by using Eqs. ([T]) and ([2]), where one replaces 


G-(k||,E)byG^B(k||,E). T 


to in-scattering processes 


lis is an approximation, since it neglects vertex corrections due 


45 


46[ |. which in general increases the conductivity. But we are 


not aware of any established method to compute those vertex corrections to linear-response 
transport within the considered framework. In our approach the Landauer transmission 
is calculated using the improved DMFT electronic structure, rather than the DFT one. 
Note that the DMFT provides the single-particle excitations of the system, whereas the 
Kohn-Sham DFT eigenvalues formally do not reveal such quasiparticle states. 

B. DMFT-based computational approach 


The transport calculations are performed according to the Green’s 
sented above by using the DFT-based transport code smeagol HQ, 


unction scheme pre- 
29l | . The many-body 
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self-energy entering in (|1]) is calculated using the EMTO-DMFT method 


34 


36 


47| 


within a screened KKR 37| approach. In both codes the Perdew-Burke-Ernzerhof (PBE) 
GGA ^ for the exchange-correlation density functional is used. Self-consistent DFT cal¬ 


culations are performed separately in smeagol and in the EMTO-code. The many-body 
self-energy is then evaluated after self-consistency in the EMTO-code, and passed to the 
SMEAGOL Green’s function to compute the transmission along Eq. ([2]). The smeagol im¬ 


ports the DFT Hamiltonian from the siesta code 


33|, which uses pseudopotentials and 


expands the wave functions of valence electrons over the basis of numerical atomic orbitals 
(NAOs). The EMTO code, in its turn, uses the muffin-tin construction; we present a de¬ 
tailed description of the projection of quantities such as the many-body self-energy from the 
EMTO basis set into the NAG basis set (smeagol/siesta) in the Appendix. 

For the EMTO-DMFT calculations, the following multi-orbital on-site in¬ 
teraction term is added to the GGA Hamiltonian in the EMTO basis: 
\Y.i{rn,a}Umm'rn''rn'''c\^a4rn'a'^irn'''a'Cirn''a- Here, Ci^a{c\ma) destroyS (creates) an elec¬ 
tron with spin a on orbital m at the site i. The Goulomb matrix elements Umm'm"m"' 
are expressed in the standard way 49(] in terms of three Kanamori parameters U, U' and 
J. Then, within DMFT the many-body system is mapped onto a multi-orbital quantum 
impurity problem, which corresponds to a set of local degrees of freedom connected to a bath 


and obeys a self-consistently condition 


23, 3- 


In the present work the impurity problem 


is solved with a spin-polarized T-matrix fluctuation exchange (SPTF) method 


26 


ml 


This method was hrst proposed by Bickers and Scalapino ^ in the context of lattice 
models. In practice, it is a perturbative expansion of the self-energy in powers of U, 
with a resummation of a specihc classes of diagrams, such as ring diagrams and ladder 
diagrams. The expansion remains reliable when the strength of interaction U is smaller 
than the bandwidth of the bath, which is fulhlled in the case of Gu-Go-Gu heterostructures 
and for the considered values of the Goulomb parameters. The impurity solver we use 
is multi-orbital, fully rotationally invariant, and moreover computationally fast, since it 
involves matrix operations like inversions and multiplications. The perturbation theory 
can be performed either self-consistently, in terms of the fully dressed Green’s function, 
or non-self-consistently, as was done in the initial implementations 1^. When the 
interaction is small with respect to the bandwidth, no appreciable difference exists between 


the non-self-consistent and self-consistent results 


53l. l54l | . This has to be distinguished from 
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the DMFT self-consistency, which is employed in both cases. In the present calcnlation, 
we nse non-dressed Green’s fnnctions to perform these inhnite snmmation of diagrams. 
Moreover, we consider a different treatment of particle-hole (PH) and particle-particle (PP 
channels. The particle-particle (PP) channel is described by the T-matrix approach ^ 
which yields renormalization of the effective interaction. This effective interaction is used 
explicitly in the particle-hole channel; details of this scheme can be found in Ref. 5^. The 
particle-particle contribution to the self-energy is combined with the Hartree-Fock and the 
second-order contributions [5^. The many-body self-energy is computed at Matsubara 


frequencies Un = (2n -|- l)n/(3, where n = 0,1, 2,... and (3 is the inverse temperature. The 


Fade 


57l | analytical continuation is employed to map the self-energies from the Matsubara 


frequencies onto real energies, as required in the transmission calculation. Note that 
since some parts of the correlation effects are already included in the GGA, the double 
counting of some terms has to be corrected. To this end, we start with the GGA electronic 
structure and replace the obtained S^^(i?) by S^^(i?) — S^^(O) in all equations of the 

], the energy E here being relative to the Fermi energy. This is a 


GGA+DMFT method 


common double-counting correction for treating metals; a more detailed description can be 


found in 


59|. 


C. Cu-Co-Cu heterostructure setup 

The basis set used in the siesta and smeagol calculations is of “double-zeta with 
polarization” (DZP) quality. In the “standard” Siesta basis construction algorithm, there 
is an “energy shift” parameter which allows to control the extent of basis functions on 
different atoms in a multi-element system in a balanced way; in our case this parameter was 
taken to be 350 meV, resulting in basis functions extending to 6.17 oq (Gu4s), 3.39 oq (CuSd), 
and 6.31 oq. (Go4s). For Co3d states, a smaller basis function localization was intentionally 
imposed, corresponding to the extention of 4.61 Oq, where Oq is the Bohr radius. The basis 
functions are usually freely chosen and not subject to optimization; however, in view of their 
quite restricted number in Siesta, it often makes sense to look at resulting ground-state 
properties of materials as a benchmark for the validity, or sufficiency, of a basis. With the 
above settings, the relaxed lattice parameters of pure constituents were found to be a = 
3.65 A (fee Gu, 1% larger than the experimental value) and a = 2.52 A, c = 4.06 A (hep 
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Co, both within 0.5% of experimental valnes). Moreover the magnetic moment per Co atom 
was 1.65 pb (equal to experiment). 

In order to calculate the transport properties, semi-inhnite leads are attached on both 
sides of the scattering region. We consider Cu(lll)-cut leads, characterized by the ABCABC 
atomic plane repetition into the leads, i.e., along the transport direction henceforth referred 
to as z. It is assumed that the ABCABC layer sequence is smoothly continued throughout 
the scattering region, including the Co monolayer (see Fig. [1]). To model the scattering 
region the sequence is repeated and the Co layer is considered to replace the Cu layer. 

In SMEAGOL the Hamiltonian of the scattering region is matched to that of the leads at 
the boundary of the scattering region, thus implying that whatever perturbation is induced 
by a scatterer it has to be conhned witin the scattering region. In other words, the simulation 
cell needs to contain enough Cu layers on each side of the Co layer to “screen” it completely. 
We verihed that using seven layers on each side provides a good agreement between the 
potential at the boundaries of our setup with the one from the periodic Cu leads calculation. 
The resulting cell geometry is shown in Fig. [1] the “Lead” + “Scattering region” + “Lead” 
composes the smeagol cell, merging on its two ends with the unperturbed semi-inhnite Cu 
electrodes. A restricted spatial relaxation within the scattering region was done by Siesta, 
whereby the total thickness of this region varied in small steps, and the z-positions of Cu2, 
Cu3 and Cu4 layers between the limiting Cul and the central Co layer were adjusted till the 
forces fell below 0.01 eV/A. This resulted in interlayer distances of 2.119 A (Cul-Cu2 and 
Cu2-Cu3), 2.118 A (Cu3-Cu4) and 2.104 A (Cu4-Co). The length for the supercell along 
the z-direction resulting from the minimization of the cell total energy was 31.619 A. The 
relaxed structural parameters obtained within the GGA were then used in the GGA-I-DMFT 
calculation, and no additional structure relaxation was attempted at the GGA-I-DMFT level. 


III. RESULTS 

In this section we discuss the changes in the electronic structure and in the conductance 
of a single Go layer sandwiched between semi-inhnite Gu electrodes caused by the inclusion 
of the Goulomb interaction at the GGA-I-DMFT level. The chosen values for Goulomb 
and exchange parameters for the 3d-Co orbitals are U = 3 eV and J = 0.9 eV, while 
no interaction beyond GGA is considered for the 3d-Cu states neither in the scattering 


DM FT supercell 



FIG. 1: (color online) Schematic representation of the supercell used in calculations. The length 
of the cell and the distances between Cul, Cu2, Cu3, Cu4 and Co are given in the text. The unit 
cell dimensions of the leads are kept at the relaxed bulk value 3.65 A. 


region nor in the leads. The values of U and J are sometimes used as htting parameters, 
although it is possible, in principle, to compute the dynamic electron-electron interaction 


matrix elements with good accuracy 


60l | . The static limit of the energy-dependent screened 


Coulomb interaction leads to a parameter in the energy range between 2 and 4 eV for all 
Sd transition metals, with substantial variations related to the choice of the local orbitals 


6l| . As the J parameter is not affected by screening it can be calculated directly within 


LSDA; it turns out to be about the same for all Sd elements, J ~ 0.9 eV The sensitivity 
of results to U and J will be briefly addressed towards the end of this section. As regards the 
temperature, two values T = 80 K and 200 K are addressed. Smaller temperatures could be 
considered; however, this would strongly increase the computational efforts connected with 
the analytical continuation of the data onto the real axis. 


A. Electronic structure calculations 

The total density of states (DOS) for the Cu-Co-Cu heterostructure is shown in Fig.[2](a). 
The Co contribution to the total DOS, attributed to the atomic sphere radius of 2.69 Oq, 
is presented in Fig. [2](b). For comparison we plot in Fig. [2](c) the Co DOS in the bulk fee 
structure, with the same atomic sphere radius. 

We start with discussing the features of the electronic structure of bulk Co. For the 
majority-spin electrons, the GGA DOS (Fig. (St) is fully occupied. In the minority spin- 
channel, the Fermi level falls between two pronounced peaks at ~ Ep±l eV. The orbital 
occupations are shown in Tab. [H According to the GGA results the majority spin-up channel 
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FIG. 2: (Color online) Densities of states calculated by EMTO in relaxed geometries. Dashed black 
lines: GGA results; solid red lines: GGA+DMFT results, (a) Total DOS per scattered region; (b) 
local DOS per central Co atom of the heterostructure; (c) local DOS per atomic sphere of the same 
size in pure fee bulk cobalt. The values of total and local Co magnetic moments for the scattered 
region are indicated in (a). 


has a nominal d-ocenpation of 4.70 while for minority electrons the ocenpation amounts to 
2.87. The s-electrons carry a negligible polarization, while p-electrons are slightly spin- 
polarized with a sign opposite to the main d-polarization which establishes a magnetic 
moment of 1.74 As a consequence of the local Coulomb interactions parameterized by 
U and J within DMFT, the DOS distribution changes considerably. The overall broadening 
is strongly modihed by the imaginary part of the self-energy. The top of the occupied d-band 
in the majority spin channel is shifted closer to the Fermi level, and some redistribution of 
the spectral weight occurs. These changes do not noticeably affect the occupation of s- 
orbitals, however, the magnetic moment, mostly due to d electrons, is signihcantly reduced 
from 1.74 /tb to 1.41 /ib. 

Essentially the correlation effects are determined not only by the magnitude of the local 
Coulomb parameters (D, J), but also by the orbital occupations. It was argued j^, 631 
that electronic interactions may lead to the creation of either a majority-spin or a minority- 
spin hole. As the majority-spin channel is essentially full, there is effectively no space for 
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TABLE I: Orbital occupations and magnetic moments for the Cu-Co-Cu heterostructure and bulk 


Co-fcc. The DMFT calculations have been performed for T = 200 K, [/ = 3 eV, J = 0.9 eV. 


^GGA 

MGGA 

^DMFT 


^^DMFT 

Atom s(4/t) P(i/t) 

da/t) (pb) Ki/t) 

Fd/t) 

dd/t) 

(he) 


Co bulk-fcc: 




Co: (0.34/0.33) (0.39/0.31) (2.87/4.70) 1.74 (0.34/0.33) (0.38/0.35) (3.04/4.50) 

1.41 


CU 4 C 0 CU 4 scattering region: 

Cul-3: (0.36/0.36) (0.33/0.33) (4.78/4.78) 0.00 (0.39/0.39) (0.45/0.45) (4.72/4.72) 0.00 

Cu4: (0.36/0.36) (0.36/0.33) (4.76/4.79) 0.00 (0.38/0.37) (0.42/0.39) (4.58/4.65) 0.02 

Co: (0.33/0.32) (0.33/0.31) (2.96/4.63) 1.63 (0.32/0.33) (0.33/0.37) (2.70/4.13) 1.48 


excitations just across the Fermi level. On the contrary, in the minority-spin channel one 
hnds a high density of electrons which can be immediately excited, leaving back holes. Such 
an occupation asymmetry has consequences concerning possible interaction channels in the 
multi-orbital Hubbard model: A majority-spin hole can only scatter with opposite-spin 
particles, which would cost an effective interaction U, while a minority-spin hole may also 
scatter with parallel-spin particles with the effective interaction U — J < U 6^. Therefore 
correlation effects are expected to manifest themselves differently for majority- and minority- 
spin electrons. 

Many DOS features of Co in the heterostructure geometry (Fig. [2 ]d) resemble those of 
bulk cobalt (Fig. [2}::), the occupation numbers of which are also given in Tab. [H As expected, 
the spin polarization in s- and p-channels is very small; moreover it is opposite to the d- 
electrons, yielding an overall magnetic moment of 1.63 fiB- As compared to the GGA case 
for the central Go layer of the heterostructure, in GGA-I-DMFT the Go s- and p-electron 
spin polarization changes sign, and the d-electrons spin splitting decreases. In the Gu-Go-Gu 
heterostructure geometry the Co-d orbitals experience hybridization with the neighboring 
Gu-d orbitals which leads to a change in the DOS of the Gu-layer in the vicinity of the 
Go-layer. 

Figure [3] depicts the local DOS of the Gu atom closest to the Go monolayer (indicated 
Gu4 in Tabled]). Even as no on-site interaction terms have been added to the 3d-Gu states, 
the Go self-energy has a large impact on the GGA-DMFT density of states. In fact, the 3d- 
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FIG. 3: (Color online) Orbital resolved DOS for Cu4 in the CU 4 C 0 CU 4 heterostructure computed 
within the GGA (black-dashed line) and GGA-I-DMFT (red solid line) for [/ = 3 eV and J = 0.9 eV. 

Cu4 states are strongly coupled with the correlated 3d-Co states and are dragged towards 
the Fermi energy, thus increasing the hybridization with the 4s- and 4p-Cu4 states. In 
contrast, the three outmost (from Co) copper layers Cu3, Cu 2 , Cul have very similar orbital 
occupations which slightly differ from those of the Cu4. The inclusion of interaction in the 
spirit of DMFT has only a slight effect upon the orbital occupations within Cu3, Cu2, and 
Cul as compared to Cu4 (see Table [Tl). Essentially the 4s- and 4p-Cu4 orbitals slightly 
increase in occupation, while 3(i-orbitals are depleted accordingly. At the same time, the 
minority/majority spin contrast gets enhanced: from 4.76/4.79 in GGA to 4.58/4.65 in 
GGA-I-DMFT. The spectral weight transfer in the Co layer - a consequence of electron 
correlations - modifies slightly, through d — d hybridization, the spin asymmetry in d-holes 
of the closest copper layer inducing a magnetic moment. 

We note that the temperature dependence of the DOS is negligible, and the effectiveness 
of electronic correlations is not signihcantly different in the heterostructure, as compared to 
the case of pure-Co fee bulk. 

These changes are typical for correlation effects in transition metals, where the self¬ 
energy near the Fermi level has Fermi-liquid character: for the imaginary part, we have 
—ImS^^Q,(i7) oc whereas the real part has negative slope, SReS^g „(i7)/cli7 < 0. Here 

E is the energy relative to the Fermi level, and a numbers the three groups of d states 
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TABLE II: Effective mass enhancements [m *for different bands of d-symmetry, calculated 
according to Eq. ([5]) as a function of the Coulomb and the exchange parameter, U and J. 


U (eV) J (eV) 

xVi 


xyy- yz^ 

1 

0.3 

1.551 

1.703 1.617 

1.578 1.737 1.656 

2 

0.6 

1.625 1.797 1.698 

1.654 1.837 1.743 

3 

0.9 

1.661 

1.791 1.703 

1.697 1.856 1.753 


in hexagonal symmetry, (z^), {xz^z) and — y‘^,xy). From the self-energy we can also 


evaluate the mass enhancement 


38j, which within DMFT amounts to 


/ jk \ f 

( m \ 


d 


I — I = 1 - S'"’" 

\^bj ^ 


QE -mb(^)) 


(5) 


where mb represents the band-mass obtained within the GGA calculations. 

The values are given in Table UTl and we note that the enhancement factors for all orbitals 
are similar, in the range of 1.6 — 1.8, which indicates that the system is medium-correlated. 


B. Transport properties 

Turning to transport properties, we display in Fig. HJ^a) the total and spin-resolved trans¬ 
mission probabilities computed with GGA and GGA-I-DMFT. The spin-resolved transmis¬ 
sion probability, T„{E), is obtained from the k-dependent transmission, Eq. ([2]), by integrat¬ 
ing over all ky-points, so that T^{E) = '^k||To-(k||, E). By inspecting Fig. HJ^a) it can 

be seen that the overall transmission is a smooth function of energy, and has a rather large 
value of about 0.5 /h in both spin channels for most considered energies, which reflects 
the fact that we deal with an all-metal junction. In GGA the transport is mainly dominated 
by the Gu-4s, -4p states, which are transmitted across the Go layer passing through the 
Go 4s states, while the Go Sd states do not contribute signihcantly to the transmission in 
this energy range. The Gu 3d states contribute to the transmission only at energies below 
— 1.5 eV. Note that the GGA-I-DMFT transmission is always smaller than the GGA one. 
The black arrows in Fig. HJ^a) indicate the energies at which a signihcant departure between 
the GGA-I-DMFT and the GGA transmission is observed. Specihcally, the spin down trans¬ 
mission drops at about 0.3 eV below the Fermi level, where the Go 3d DOS is high in the 
DMFT results, while the GGA transmission stays rather constant. In contrast, the spin-up 
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transmission shows a “bnmp”, which extends over a region of abont 2 eV aronnd the Fermi 
level. The slight dip within this bnmp at abont —0.5 eV is at the same energy as the peak 
in the Co 3d DOS, and represents a Fano-type rednction of transmission in snch a metallic 
system dne to interference of electrons in different condncting channels {m, 3, 651. In onr 
calcnlations this featnre is a conseqnence of electronic correlations on the Co atom, which 
throngh the d — d hybridization indnce spin-polarization of the Sd-Cnd states and simnlta- 
neonsly prodnce the shift in the DOS of Fig. [3l In general, we note that for snch all-metal 
systems the relation between DOS and transmission is non-trivial, since interference effects 
can lead to enhanced transmission also for energies with low DOS; alternatively, for high 
DOS the increased nnmber of pathways for electrons can lead to a decrease of transmission. 

From a many-body perspective, the added self-energy contribntes in dephasing the elec¬ 
trons dnring the flow throngh the scattering region, so that the Landaner transmission 
compnted with the many-body Green’s fnnction is expected to be reduced in comparison 
with the DFT case. In principle, the opposite effect, namely the effective in-flow of elec¬ 
trons from the many-body self-energy “electrode” into the scattering region, would tend to 
increase the transmission. However, this in-flow process is not included in our calculations. 


as it is related to the vertex corrections 


45|. 


The spin polarization of the transmission is computed according to the formula 


( 6 ) 


T,{Ep)+Ti{EFy 

for either DFT, or DFT-I-DMFT, where is the transmission for the spin channel a. 

The spin polarization in transmission obtained by GGA yields = 0.18 [see also 

Fig. Ill[b)], while the GGA-I-DMFT value reaches 0.33 at the Fermi level, and increases up to 
almost 0.8 at slightly lower energies. These results demostrate that electronic correlations 
may be decisive and lead to an increase in the spin polarization of transmission. As seen 
in Fig. Il][b), the enhancement in the spin polarization with respect to the GGA result 
is essentially temperature independent in the energy range of Ep ± 0.5eV. Therefore we 
conclude that the enhanced spin contrast in transmission is a many-body effect rather than 
a temperature fluctuation effect. 

Finally, we test how the results for the transmission depend on the strength of the local 


Goulomb interaction parameters U and J. The interaction matrix e 
usually parametrized using Slater integrals (F^) with k = 0,2,4 


ements Umm'm"m"' sre 


49l |. Accordingly the 
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FIG. 4: (Color online) Left: Spin-resolved transmission: majority-spins F]- (upper panel), mi¬ 
nority spins (middle panel) and total (lower panel). Black-dashed/red-solid lines are the 
GGA/GGA-I-DMFT results. Right: Transmission spin polarization around the Fermi energy GGA 
(black-dashed) and at finite temperatures T=80K (red-dotted-dashed), T = 200 K(blue solid). 
The Coulomb and exchange parameters are U = 3 eV, J = 0.9 eV. 


Hubbard U parameter is constructed as a simple average over all possible pairs of correlated 
orbitals and is identified with the Slater integral U = F^. The other Slater integrals 
are htted to the multiplet structure measured in X-ray photoemision 6^. An empiric 


relation has been introduced whic 
Slater integrals, F^/F‘^ ^ 0.625 


1 connects the magnitude of the second and forth order 


67| . The Hund’s exchange J is expressed in terms of and 


F^ which for the d-shell takes the form J = (F^-|-F^)/14 68|, therefore the knowledge of the 
(t/, J)-pair allows to compute all the matrix elements Umm'm"m"' ■ In Fig. [5](a) we plot the 
transmission (computed at 200 K) keeping the ratio U/J = 1.0/0.3 constant. While scaling 
the ratio a ■ U/J with a = 1,2,3 we observe a monothonic reduction of the transmission 
at the Fermi level. This result is expected, as the matrix elements of the interaction are 
scalled in magnitude. In the same time, larger mass enhancement factors are obtained as 
a increases (see Tab. [I]). Consequently we may conclude that the heavier the electron is, 
the smaller is the transmission at the Fermi level. Within ±0.15 eV of the Fermi level, a 
flat region in the transmission can be seen. Beyond these ranges, we note that below the 
Fermi level, down to —1 eV from it, the transmission decreases almost indiscriminately for 
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FIG. 5: (Color online) Transmission for different values of Coulomb parameters: (a) with a fixed 
ratio U/J, (b) fixed J and (c) fixed U, ai T = 200 K. 

different values of {U, J). In the positive energy range up to roughly +1 eV, on the contrary, 
the transmission values differ, with larger {U, J) resulting in lower transmission. 

In Fig. [5](b) we display the dependence on parameters differently keeping J = 0.9, eV 
constant and varying U. We note that with higher U, the flat region centered at Ep shrinks 
a bit, which can be traced back to a stronger presence of d orbitals in the correlated trans¬ 
mission. The inset of Fig. [5](b) depicts the “contrast” or spin difference in transmission: 
AT{E) = T-^{E) — T^{E). The spin contrast changes slope as U varies, from markedly 
positive d{AT)/dE at Ep for f/ = 2 eV to roughly zero for f/ = 3 eV. 

In Fig. [5](c), on the contrary, the U parameter is kept hxed to a “good” (yielding large 
flat region) value U = 3 eV, and two different values are taken for the J parameter. On 
increasing J, the flat region around Ep gets further reduced. Simultaneously (as seen in the 
inset and contrary to the behavior depicted in Fig. [5 )d), AT{E) acquires a negative slope. 
The change in the slope is potentially an interesting effect to be exploited in thermoelectric 
transport. 
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IV. CONCLUSION 


In this work we studied the effects of local electronic interactions and finite temperatures 
upon the transmission across the CU 4 C 0 CU 4 metallic heterostructure. Electronic structure 
calculations were performed using GGA and GGA+DMFT, assuming a local Goulomb in¬ 
teraction in the Go layer. We used a fully rotationally invariant Goulomb interaction on 
cobalt d-orbitals. The effective mass enhancement ratio for all orbitals is in the range of 1.6 
to 1.8, suggesting that we deal with medium correlated system. Goncerning the density of 
states, the presence of Goulomb interactions leads to a shift of the majority-spin channel of 
Go d-orbitals towards the Fermi level and to a redistribution of the spectral weights. In the 
minority spin-channel, the changes are less pronounced. This difference leads to different 
correlation effects for the majority- and minority-spin electrons. All these causes a decrease 
of the overall Go magnetic moment (with a predominant d character) from 1.63 (GGA)to 
1.48 fiB (GGA+DMFT). 

The transmission probability has been computed combining two different ab-initio codes, 
EMTO and Siesta/Smeagol. In order to transfer the many-body self-energy computed 
within the EMTO code into Siesta/Smeagol, we used a nearly unitary transformation 
which can be determined by requiring that the expectation value of the occupation matrix 
should be representation independent. The methodology was illustrated for a two-orbitals 
model (see Appendix [B]) and was carried out numerically for the GU 4 G 0 GU 4 heterostructure. 
Note that such a transformation is rather general and can be used in transfering quantities 
between two different implementations. Several tests confirmed that the proposed method is 
robust and numerically stable. With this combination of methods we have studied the trans¬ 
mission as a function of temperature and Goulomb parameters, which reveals the metallic 
character of the system considered. Substantial differences in the conducting processes, re¬ 
lated to the presence of local Goulomb interaction in the cobalt layer were observed, due to 
changes in the electronic structure. Generally the transmission decreases with interaction, 
although the relation between changes in the electronic configuration and the transmission 
is highly non-trivial, due to interference effects. With electron correlations properly taken 
into account, the total transmission at the Fermi level drops by about 20%, whereas its spin 
polarization (spin contrast) increases by about 40%. These effects are entirely a consequence 
of the electronic correlation, since the transmission is practically temperature independent 
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in the range of Ep dz 0.5 eV. This suggest that the enhanced spin contrast in transmission is 
fully a many-body effect. In order to quantify the spin polarization effects in the transmis¬ 
sion, we studied the transmission difference AT{E). This quantity clearly displays a strong 
dependence on the Coulomb parameters. 

In conclusion, we have shown that electronic correlations may considerably affect the 
transmission and spin hlter properties of heterostructures, even though the correlations 
would be classihed only as “medium” when considering the effective mass enhancement. 
Hence results based on studies neglecting electronic correlations, which are numerous, must 
be interpreted with caution. 
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APPENDIX 

Appendix A: SIESTA and EMTO density of states 

Inorder to demonstrate the reliability of both codes concerning the electronic structures, 
we present below the density of states for the heterostructure (Fig. E]). A rather good 
agreement is apparent. 
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E-Ep(eV) 


FIG. 6: (color online) Total (unit cell) and Co density of states computed within the Siesta and 
EMTO technique. 

Appendix B: Matrix elements of the self-energy in the NAO basis set 


The matrix representation of the self-energy operator is determined by the chosen ba¬ 
sis. Within the EMTO basis set it has the form ^j;^,( 2 ;). Here R and B! are site 

indices while the L symbol labels the orbital quantum numbers. To compute the transmis¬ 
sion/conductance it is desirable to work within the Siesta/Smeagol basis. We emphasize 
that the major part of calculation is done within the Siesta+Smeagol package. 

Since signihcant methodological differences exists between Siesta and EMTO, in the 
following we discuss the methodology of data transfer between the two codes. We describe 
briefly the most signihcant differences. The former (Siesta) implementation uses norm- 
conserving pseudopotentials, whereas the latter (EMTO) uses an all-electron formulation. 
Siesta uses no shape approximation with respect to the potential, whereas EMTO relies 
on the muffin-tin concept 3J, [3^. Basis functions in Siesta are atom-centered numerical 
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functions, whose angular parts are spherical harmonics, and radial parts are strictly confined 
numerical functions. A tradition hnding its origins in quantum chemistry suggests that, in 
order to improve variational freedom of the basis set, more than one radial function is 
adopted into the basis for a given angular combination {l,m), referred to as “multiple-^” 
basis orbitals. 

Even as Siesta maintains flexibility in constructing the basis set out of different “zetas”, 
possibly including moreover “polarization orbitals” and allowing a variety of schemes to 
enforce conhnement, we shall stick in the following to the case of just “double-^”, i.e., 
2x5 = 10 basis functions, provided to accurately describe the 3d states on each cobalt atom. 
We shall fix some notation for further reference. The cumulate index of a basis function 
will be fi={nm(}, where I indicates the atom carrying the basis function, ( numbers the 
“zeta”s (= 1 or 2 in our case), and the {l,m) is the conventional angular moment index. It 
should be noted, however, that Siesta employs real combinations of “standard” spherical 
harmonics, so that the indices m = —2 through 2 for I = 2 correspond to xy, yz, —r^, 
xz and x'^ — y'^ d-functions, correspondingly. With the above notation, the Ah eigenstate of 
the Kohn-Sham Hamiltonian will expand into the basis functions 0^ as follows: 


-Rj), 


(Bl) 


and the variational principle yields the expansion coefficients SSj. Within the EMTO, the 
d-orbitals manifold is constructed using a basis set with real harmonics (physical orbitals 
representation), and the occupation matrix is obtained integrating the complex contour 
Green’s function (properly normalized path operator) in terms of the exact muffin-tin orbitals 


r) corresponding to the energy e 34j-l36l|. On the other hand, in Siesta the numerical 


basis set is not restricted to the physical orbitals and allows the definition of simple/double 
( “atomic-like” representations. Enhancing the numerical atomic orbitals basis vectors does 
not affect the dimension of the vector spaces (in both cases the d-subspace), however, it 
complicates the algebra for the transformation matrix. The transformation involves the 
double zeta basis set, 0^, which for the L(/, m)-manifold contains 2(2/ -|- 1) components 
splitted into the first (2/ -|- 1)—single-zeta and the second (2/ -|- 1)—double-zeta components. 
First we write the explicit form for the EMTO orbital in a vector form for the L{1 = 2,m = 
2/-I-1) subspace |TijL(e, r)), and the corresponding Siesta basis vector in the double C-basis, 
|0HL(r)) = (0^‘’(r)) 0^'’(r))- The general transformation matrix is defined as the 
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inner products of \I/H'L'(e, r) with Siesta’s double-^ basis, and takes the form of a dyadic 
product: 


_ J 


(B2) 


The dehnition Eq. flB2p for the V matrix suggests the possibility of an explicit construc¬ 
tion. However, a couple of obvious ambiguities may arise: (i) the transformation carries an 


energy dependence originating from the energy dependence of the EMTO orbitals 


34 


in contrast to the NAO basis set that is energy independent; (ii) normalization of the scat 


tering path operator of EMTO is performed in a particular screened representation 


34 


a more involved procedure in comparison with the straight normalization of the NAO basis 
set. 

Moreover, it is important to keep in mind that the closure relations that are typically 
used to build the matrix transformations are valid on the Hilbert space spanned by the 
eigenvectors of the Hamiltonian. These relations hold exactly for the numerical results of 
each code separately; nevertheless, the numerical results produced by different codes are not 
identical in the mathematical sense. Small differences between the observables computed 
with different codes may exist due to various factors - from unwanted numerical roundoff 
errors to incompleteness of the basis sets. This indicates that the closure relation for the 
Hilbert space of the EMTO will not be exact (in the mathematical since) when used in 
the Hilbert space spanned by the eigenvectors of Siesta. Consequently, we expect that 
the transformation matrix will fulhll the usual requirements (i.e., unitarity) only within 
numerical inaccuracies. 

In view of these formal difficulties in obtaining a basis transformation to match a multiple 
scattering method with a Hamiltonian based scheme we propose an approach using the fact 
that expectation values should be independent of the specihc representation. We apply 
this fundamental concept of quantum mechanics to the orbital occupation matrix (density 
matrix, for the d-manifold (i.e., rrii, rrij = 1,..., 5), and we look for a formally similar 

matrix transformation W : 


rrij ,mi mj ,mi J 


EMTO ( mi,rrij | ^ SIESTA 

'‘'rriimj \ ) '^mi,mj 


(B3) 
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Equation (lB3j) provides us with a system of equations to determine the matrix elements 
of W numerically. Indeed, by inspecting the diagonal elements of the occupation matrix 
for the xy, yz, — r^, xz and x"^ — y"^ orbitals respectively, we obtained the data shown 
in Table Hm One should note that in the physical-EMTO basis the occupation matrix is 
diagonal. The corresponding Siesta occupation matrix has non-diagonal elements that are 
about two orders of magnitude smaller than the diagonal ones. While the symmetry and the 



rrii = xy 

yz 

z2 _ J.2 

zx 

to 

1 

to 

EMTO 

t 

0.629 

0.526 

0.658 

0.526 

0.629 

i 

0.922 

0.925 

0.908 

0.925 

0.922 

^Siesta 

^TTii ylTLi 

t 

0.639 

0.485 

0.784 

0.485 

0.639 

i 

0.935 

0.952 

0.924 

0.952 

0.935 


TABLE III: Occupation matrix of Co-d orbitals in the Co/Cu heterostructure 

qualitative trends in the occupations are the same, the exact numerical values are not. In 
other words, Eq. flB3l) is an approximative (numerical) relation, however, the resulting W 
matrix reflects the symmetry of The non-zero elements on the diagonal are close to 

1 for the £rst-C block and take small imaginary values for the second-^ one. The equations 
Eq. fIClj) and flC4D provide explicit values. 

Accordingly, given the matrix elements for the transformation matrix, the self-energy 
generated in the EMTO-basis set can be transferred to the double(multiple)-zeta Siesta 
basis set according to 

^SIESTA ^ (-34) 

In the following, we give simple examples to illustrate the self-energy transformation from 
the physical into a numerical double-zeta basis. 


1. Example: a two orbital model in the cubic symmetry 


We consider a simplihed case of a diagonal self-energy corresponding to a two orbital 
model in the cubic symmetry. The self-energy and the “occupation” matrix can be written 
as: 


Y^emto 





. EMTO 

. f i 



(B5) 


22 




















In such a case all inner products of the type r)|0m^ (r)) are zero unless m* = rrij, 

such that the transformation matrix has a generic form: 


mi,mi \ / 


f 0 \ 

0 hh22^ 


w, 


2( 


11 


0 


0 IlfP 


Re-writing Eq. (163^ : 
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V 0 
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0 
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one hnds 
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0 




V 0 


\Ai2 
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and accordingly the self-energy in the NAO basis-set has the form: 

/ 


E = 


El 


^ ni 

0 

^ ni 

0 ^ 

0 

^2 — 

^ n2 

0 

2 712 

vww 

ni 

0 

- 2 C 

^ ni 

0 

0 

/- 1 C- 2 C 

2 n2 

0 

-2C 

y' ^2 , 

^2^ / 


/ SN’R 

S1C,2C ^ 


S2C,2C j 


(B6) 


(B7) 


(B8) 


(B9) 


in which every is a 2 x 2 block-diagonal matrix. There are a couple of conclusions 

to be drawn from the above simplihed example: (i) as a consequence of the reduced weight 
in the second ( occupation (at least 10 times smaller that in the single C), the magnitude 
of the matrix elements of the self-energy in the NAO basis set follow the relation > 

^1C,2C ^ ^2C,2C. 

(ii) the existence of only-diagonal orbital occupations does not imply the 
existence of an unitary transformation, except for the case of numerical identical values 
for the occupation matrices in both basis. These observations do not change when the 
symmetry is lower than cubic and non-zero matrix elements on the non-diagonal of the 
occupation matrix appear. 
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Appendix C: Transformation matrix for the Co-d manifold in the Cu-Co-Cu het- 


erostructnre 


This section provides the spin-resolved matrix elements, used in the calculation. 


/ 1.024 0 0 

0 0.987 0 

0 0 1.096 


W-f = 


V 


0 

0 

0.178i 
0 
0 
0 
0 


0 
0 
0 

0.218i 

0 

0 

0 


0 

0 

0 

0 

0 

0 

0 


0 

0 

0 

0.985 

0 

0 

0 

0 

0.218i 

0 


0 ^ 
0 
0 
0 

1.024 
0 
0 
0 
0 

0.178i y 


The orthogonality can be checked using the following relations: 

/ 1.016 0 0 0 0 \ 

0 0.922 0 0 0 


= 


and 


ITtWy = 


/ 1.048 
0 
0 
0 
0 

0.182i 
0 
0 
0 
0 


0 

0.970 

0 

0 

0 

0 

0.215i 

0 

0 

0 


0 

0 

1.202 

0 

0 

0 

0 

0 

0 

0 


0 

0 

0 

0 

0 

0 

0.970 

0 

0 

0 

0 

0.215i 

0 


0 
0 
0 

0 

0 

0 

0 

1.047 


1.202 0 0 
0 0.922 0 

0 0 1.016 / 


0.182i 

0 

0 

0 

0 


0 -0.032 

0 0 

0 0 

0 0 

0.183i 0 




0 0 0 0 

0.215i 0 0 0 

0 0 0 0 

0 0 0.215i 0 

000 0.182i 

0 0 0 0 

-0.048 0 0 0 

0 0 0 0 

0 0 -0.048 0 

0 0 0 -0.032 J 


(Cl) 


(C2) 


(C3) 
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The corresponding transformation matrix for spin-down component reads: 


in 


^ 0.960 0 0 -0.007 0 ^ 

0 0.9572 0 0 -0.007 

0 0 0.967 0 0 

-0.007 0 0 0.957 0 

0 -0.007 0 0 0.960 

0.304 0 0 -0.002 0 

0 0.337 0 0 -0.002 

0 0 0.287 0 0 

-0.002 0 0 0.337 0 

^ 0 -0.002 0 0 0.304 y 




/ 1.014 0 0 -0.015 0 \ 

0 1.029 0 0 -0.015 

0 0 1.017 0 0 

-0.015 0 0 1.029 0 

\ 0 -0.015 0 0 1.014 / 


(C4) 


(C5) 




/ 0.922 0 0 -0.014 0 0.292 0 0 0 0 ^ 

0 0.916 0 0 -0.014 0 0.322 0 0 0 

0 0 0.935 0 0 0 0 0.278 0 0 

-0.014 0 0 0.916 0 0 0 0 0.322 0 

0 -0.014 0 0 0.922 0 0 0 0 0.292 

0.297 0 0 0 0 0.092 0 0 0 0 

0 0.325 0 0 0 0 0.114 0 0 0 

0 0 0.278 0 0 0 0 0.083 0 0 

0 0 0 0.325 0 0 0 0 0.114 0 

^ 0 0 0 0 0.292 0 0 0 0 0.092 y 

(C6) 
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Appendix D: Assessment of accuracy 


In order to test the effect of the numerical inaccuracies occurring in the transformation 
matrix on the final results, we compute the total transmission by using a simplified model 
for the transformation matrix: 


Model — 




(Dl) 


where = al and = (31, where I is the unity matrix. The results for a = 0.9, 
(3 = 0.1; a = 0.5, (3 = 0.5; and a = 0.1, (3 = 0.9, respectively are given in Fig. [71 It can be 
clearly seen that even for such a crude approximation, the results for the first model (i.e., 
with a significant weight of the self-energy on the first zeta orbital) differ with only a few 
percent over large energy domains. For the occupied states, large values for provide 
already a good approximation for the transformation. 



FIG. 7: (color online) Comparison between the transmission functions obtained for different model- 
forms of IFF’2C, 
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